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Under  the  sponsorship  of  The  Office  of  Naval  Research  (Code  12)  AHR 
Program,  The  Applied  Research  Laboratory  of  Penn  State  University 
performs  basic  research  in  hydrodynamics  and  hydrodynamic  noise.  The 
hydrodynamics  research  conducted  under  this  program  falls  into  two  basic 
thrust  areas: 


To  develop  an  improved  understanding  of  the  fluid  mechanics  and 
acoustics  associated  with  low- speed  turbomachines  and  marine 
propulsors .  To  employ  this  knowledge  to  the  development  of 
improved  propulsor  and  turbomachine  design  methods.  /* 


To  develop  fundamental  understanding  of  the  mechanisms  that 
cause  drag  on  bodies  and  surfaces  and  to  explore  novel  methods 
to  reduce  drag.  .  1  ^  ^ _ 

Under  each  thrust  area,  one  or  more  projects  are  conducted  under 
the  direction  of  the  principal  investigator  who  initiated  the  given 
task.  All  tasks  are  designed  to  provide  results  that  will  improve  the 
scientific  understanding  of  various  hydrodynamic  phenomena  associated 
with  the  operation  of  submerged  bodies  and  surfaces. 


This  report  documents  the  technical  progress  realized  during  the 
first  quarter  of  FY  90  for  the  projects  currently  approved  under  this 
program. 
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TURBOMACHINERY 


SUBTASK  T1 


Title :  Turbomachine  Internal  Flow  Definition 

(S.  A.  Abdallah,  Cincinnati,  and  M.  L.  Billet.  Penn  State) 

BACKGROUND 

The  internal  flow  field  of  a  wake  adapted  turbomachine  is  dominated 
by  three  dimensional  and  unsteady  effects.  The  three  dimensionality  of 
the  flow  field  is  demonstrated  by  the  strong  secondary  flows  which  have 
been  experimentally  measured.  The  unsteadiness  of  the  flow  is  due  to 
the  interaction  of  the  downstream  blade  rows  with  the  wakes  shed  from 
the  upstream  blade  rows.  The  development  of  computational  tools  to 
predict  accurately  these  types  of  flow  fields  is  essential  if  successful 
development  of  high  performance  turbomachinery  is  to  be  achieved. 


PROGRESS 

In  order  to  participate  in  the  SSPA-CHT- 11HR  workshop  on  Ship 
Viscous  Flow,  Sweden,  12-14  September  1990,  we  changed  our  research 
goals  from  unsteady  internal  flow  to  steady  turbulent  flow  over  a  whole 
ship.  The  change  was  approved  by  Mr.  J.  Fein. 

Preliminary  results  have  been  obtained  for  three-dimensional 
turbulent  flow  over  a  wiggly  hull.  The  Baldwin  Lomax  turbulence  model 
is  used  in  the  calculations.  Comparisons  with  experimental  data  are  in 
good  agreement. 

Enclosed  is  a  paper  entitled  "The  Discrete  Continuity  i:  imitive 

Variable  Solutions  of  Incompressible  Flow."  Submitted  to  the  Journal  of 
Computational  Physics  and  published  in  the  proceedings  of  the  "Third 
International  Congress  of  Fluid  Mechanics,"  Cairo,  Egypt  (Jan.  2-4, 

1990) . 
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DRAG  REDUCTION 


SUBTASK  DR1 


Title :  Turbulent  Spot  Generation  by  Freely  Suspended  Particles 

(H.  L.  Petrie  and  P.  J.  Morris) 

BACKGROUND 

Laminar  flow  control  (LFC)  has  been  an  attractive  technology  for 
many  years  but  premature  transition,  induced  by  particles  in  the  fluid, 
has  prevented  its  successful  hydrodynamic  implementation.  Although  an 
extensive  amount  of  experimental  work  and  analysis  has  examined  the 
transition  induced  by  various  types  of  roughness  elements  and 
intermittent  disturbances  such  as  sparks,  basic  questions  regarding  the 
mechanisms  for  transition  induction  by  freely  suspended  particles  still 
remain  unanswered.  The  subject  experimental  and  analytic  study  of 
particle  induced  transition  is  concerned  with  these  issues  and  is 
entering  the  ninth  month  of  activity. 

Heated  body  LFC  experiments  have  demonstrated  clearly  the 
degradation  of  LFC  performance  due  to  freestream  particles.  However,  a 
detailed  explanation  of  this  phenomenon  has  not  been  developed. 
Experiments  by  Lauchle ,  et  al.  [1],  observed  that  nearly  neutrally 
buoyant  particles  of  sufficient  size  caused  turbulent  spots  on  an  LFC 
heated  body  when  seeded  into  the  freestream.  However,  these  particles 
failed  to  generated  spots  when  injected  through  a  dye  port  at  the  nose 
of  the  body.  If  these  particles  did  not  migrate  away  from  the  wall 
significantly  while  accelerating  around  the  nose,  then  freestream  seeded 
particles  must  act  on  the  boundary  layer  with  a  different  mechanism  than 
the  wall  injected  particles.  This  suggests  that  the  dynamics  of 
crossing  streamlines  and  possibly  wall  impact  may  be  crucial  to  spot 
formation.  In  support  of  this,  Hall  [2],  observed  that  transition  could 
be  induced  when  a  spherical  particle  supported  by  wire  was  pulled  into 
the  wall  even  though,  once  fixed  to  the  wall,  the  particle  was  too  small 
to  induce  a  turbulent  spot. 


PROGRESS 

As  reported  earlier,  the  analytical  approach  is  addressing  two  main 
issues  that  support  the  experimental  measurements.  The  first  is  the 
determination  of  particle  trajectories.  The  second  is  the  calculation 
of  th~  development  of  wave  packets  generated  by  the  moving  particles. 
During  this  reporting  period  we  have  concentrated  on  two  aspects  of  the 
analysis  related  to  the  trajectory  calculations. 

The  basic  flow  field  has  been  defined  using  boundary  layer 
analysis.  The  Blasius,  flat  plate  boundary  layer  solution  has  been 
matched  to  second  order  with  the  external  uniform  flow.  This  solution 
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includes  the  external  flow  due  to  the  displacement  thickness  of  the 
boundary  layer.  It  also  provides  a  smooth  transition  for  the  velocity 
normal  to  the  wall  from  its  non- zero  value  at  the  edge  of  the  boundary 
layer  to  its  zero  value  in  the  freestream.  Calculations  have  been 
performed  for  the  test  conditions  to  be  used  in  the  water  channel.  The 
results  will  be  compared  with  the  experimental  observations  during  the 
calibration  phase  of  the  flat  plate  model. 

With  the  basis  defined  it  is  possible  to  determine  particle 
trajectories.  The  trajectories  depend  on  the  forces  acting  on  the 
particle.  Outside  the  boundary  layer  these  forces  are  well  understood 
since  the  flow  is  uniform.  Particle  trajectory  calculations  are  being 
performed  in  this  region.  These  will  be  used  to  determine  the 
appropriate  release  location  for  particles  in  this  experiments.  Once 
the  particle  enters  the  boundary  layer  the  situation  is  more  complicated 
as  the  basic  flow  is  sheared.  In  addition  to  the  forces  present  in  the 
freestream  there  is  a  lift  force  due  to  the  shear  and  to  particle 
rotation.  The  proximity  of  the  wall  to  the  particle  can  also  influence 
the  lift.  An  asymptotic  analysis  for  low  particle  Reynolds  numbers  has 
been  performed  by  Saffman.  In  fact  the  lift  force  has  become  known  as 
the  "Saffman  lift  force."  However,  we  envision  particle  Reynolds 
numbers  of  the  order  of  300.  Several  experiments  have  been  performed  to 
determine  the  lift  force  as  a  function  of  local  shear  and  distance  from 
the  wall.  Some  experiments  involve  the  direct  measurement  of  lift, 
which  is  extremely  small.  Others  involve  a  balancing  of  the  lift  force 
with  buoyancy  force  in  an  inclined  channel.  The  uncertainty  in  the 
experimental  data  is  quite  large.  At  present  we  are  correlating  the 
data  to  develop  an  empirical  extension  of  the  Saffman  formula  to  higher 
Reynolds  number.  This  is  in  the  same  spirit  as  the  extension  of  the 
Stokes  drag  coefficient.  The  separate  effects  of  wall  proximity  will 
also  be  included.  This  empirical  model  will  be  incorporated  into  the 
trajectory  prediction  scheme  and  predictions  will  be  made  for  the 
particle  motion  within  the  boundary  layer.  These  predictions  will  be 
compared  with  the  particle  trajectory  measurements  when  they  are 
available,  and  modifications,  if  necessary,  will  be  made  to  the 
empirical  lift  force  correlations. 

Experiments  are  to  be  conducted  in  the  low  turbulence  open  channel 
flow  loop  facility  of  the  Aerospace  Engineering  Department  at  Penn 
State.  Nearly  all  of  the  various  pieces  for  the  channel  and  test  plate 
have  been  completed.  The  list  of  components  that  have  been  fabricated 
or  modified  includes  the  test  plate  with  various  subassemblies,  the 
channel  floor  and  sidewall  support  frame,  the  test  section  support 
frame,  the  diffuser,  the  settling  section  inlet,  and  the  flow 
conditioning  section  in  the  settling  chamber.  This  represents  a 
substantial  design  and  fabrication  effort  that  has  occupied  the  time  of 
numerous  personnel  and  is  a  significant  milestone  now  that  we  are 
beginning  our  ninth  month  of  project  activities.  Channel  assembly  and 
test  plate  instal 1 ation  will  commence  early  in  the  next  quarter. 

Initial  familiarization  with  the  computer  based  video  imaging 
system  for  particle  pathline  measurement  has  begun.  To  date,  this  has 


6 


been  limited  to  ensuring  that  all  components  function  as  they  should  and 
evaluation  of  the  systems  software. 

PLANS 

A  major  activity  of  the  next  quarter  will  be  the  installation  of 
the  plate,  the  assembly  of  the  channel  and  then  flowfield  evaluation. 
Initial  experimental  work  will  begin  late  in  the  quarter  with  the 
emphasis  being  on  channel  flow  quality  evaluation.  At  first,  flow 
visualization  will  be  used  to  examine  the  flow  from  the  inlet  of  the 
settling  chamber  to  its  exit  at  the  test  section.  The  goal  of  this  work 
is  to  determine  a  configuration  of  the  flow  conditioning  section  that  is 
sufficient  to  quell  the  large  disturbances  upstream  of  it.  Following 
this,  the  test  plate  sidewall  suction  slots  and  tail  must  be  adjusted. 
This  requires  both  flow  visualization  and  boundary  layer  surveys. 
Coincident  with  these  efforts,  imaging  techniques,  particle  release 
trials,  and  development  of  general  methodology  for  the  particle 
experiments  will  begin.  Work  on  the  particle  trajectory  analysis  will 
continue  with  the  calculations  being  used  as  a  guide  for  selection  of 
particle  release  positions. 

REFERENCES 

1.  Lauchle ,  G.  C.,  Petrie,  H.  L. ,  and  Stinebring,  D.  R. ,  "Effects  of 
Particles  on  the  Delayed  Transition  of  a  Heated  Body,"  ARL 
Technical  Memorandum  86-213. 

2.  Hall,  W.  R. ,  "Interaction  of  the  Wake  from  Bluff  Bodies  with  an 
Initially  Laminar  Boundary  Layer,"  AIAA  Journal.  V5 ,  N8 ,  pp.  1386- 
1392.  "Boundary  Layer  on  a  Flat  Plate,"  Proc .  Rov . . 
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SUBTASK  DR2 


Title :  Microbubble  Injection  in  Ax i symmetric  Flows  (S.  Deutsch) 


BACKGROUND 

It  has  been  demonstrated  that  the  injection  of  gas  to  form 
microbubbles  in  a  liquid  turbulent  boundary  layer  on  a  flat  plate  at 
nominally  zero  pressure  gradient  can  reduce  skin  friction  drag  by  as 
much  as  90%  locally.  The  extension  of  these  results  to  submerged 
axisymmetric  bodies  is  the  subject  of  the  task. 


PROGRESS 

A  publication  based  on  an  extensive  set  of  local  skin  friction 
measurements  has  been  accepted  (with  revision)  by  the  Physics  of  Fluids. 
Major  results  covered  in  the  paper  include  the  fact  that  the  persistence 
of  the  phenomenon  on  an  axisymmetric  body  is  the  same  as  that  on  a  flat 
plate  and  that  the  poor  performance  of  microbubble  drag  reduction  at  low 
speeds  is  a  result  of  a  buoyancy  driven  instability  that  removes  bubbles 
from  the  boundary  layer. 

We  are  continuing  our  study  of  the  effect  of  pressure  gradient.  We 
have  repeated  our  integrated  skin  friction  measurements  and  LDV  studies 
and  have  measured  the  pressure  gradient  on  the  body  with  and  without 
bubble  injection.  Although  data  analysis  is  continuing,  it  appears  that 
the  earlier  drag  reduction  data  is  repeatable  and  that  injection  of  gas 
does  not  change  the  pressure  gradient  on  the  body.  Relative  levels  of 
drag  reduction  would  appear  then  to  be  correct,  independent  of  the 
effect  of  pressure  gradient  on  the  balance  results.  Our  earlier 
conclusions  that  an  unfavorable  pressure  gradient  led  to  higher  levels 
of  drag  reduction  and  the  prospect  of  early  separation,  while  a 
favorable  gradient  led  to  much  lower  levels  of  drag  reduction  would 
appear  to  be  confirmed. 


PLANS 


Analysis  of  the  LDV  data  will  continue.  This  should  help  determine 
the  influence  of  the  pressure  gradients  on  the  drag  balance  results.  In 
addition,  we  have  begun  simple  computations  of  the  flow  field  over  the 
axisymmetric  body  for  the  zero,  adverse  and  favorable  pressure  gradient 
cases.  These  computations,  which  employ  the  measured  pressure 
gradients,  will  give  us  yet  another  test  of  the  reliability  of  the 
(absolute)  values  of  the  balance  results  in  a  pressure  gradient. 

The  results  and  plans  sketched  above  will  provide  the  basis  for  a 
M.  S.  Thesis  in  Aerospace  Engineering  for  Mr.  H.  Clark  as  well  as  for  an 
additional  publication. 
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THE  DISCRETE  CONTINUITY  IN  PRIMITIVE  VARIABLE  SOLUTIONS 
OF  INCOMPRESSIBLE  FLOW 

F.  Sotiropoules  and  S.  Abdallah 
Department  of  Aerospace  Engineering  and  Engineering  Mechanics 
University  of  Cincinnati,  Cincinnati,  Ohio  U.S.A. 


ABSTRACT 

The  use  of  a  non -staggered  computatiennal  grid  for  the  numerical  solutions 
of  the  ineenpre ss I fc 3 e  f  1  v-  equations  has  many  advantages  over  the  use  of  a 
stnr.qe’c <3  grid.  A  penalty,  however,  is  inherent  in  the  finite-difference 
.-:rri  oxirrsti  ons  of  the  governing  equations  on  ncr.- staggered  grids.  In  the 
primitive-variable  solutions,  the  penalty  is  that  the  discrete  continuity 
equation  dees  not  converge  tc  a  machine  accuracy.  It  rather  converges  to  a 
source  term  which  is  proportional  to  the  fourth  order  derivative  of  thr 
pressure,  the  time  increment  and  the  square  of  the.  grid  spacing.  \An  approach 
which  minimizes  the  error  in  the  discrete  continuity  equation  i.«  developed.. 
Numerical  results  obtained  for  the  driven  cavity  problem  confirm  the  analytical 
developments  -*  •• 
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higher  order  terms’  In  equations  (18a)  and  (16b)  respectively 
variables  defined  in  equations  (16a)  and  (16b)  respectively 
pressure  divided  by  density 
Reynolds  number 
time  • 

velocity  components  in  x-  and  y-dircctions  respectively 

spece  coordinates 

source  term  in  equation  (27) 

variables  defined  in  equations  (Aa)  ar.d  (Ab)  respectively 
time  and  space  increments 

second  order  accurate  finite-differer.ee  approximations  for  the 

partial  derivatives  3/3:<  and  3/cy  respectively 
second  order  accurate  finite- difference  approximations  for  the 

partial  derivatives  S2/cx2  and  <?2/<?y2  respectively. 


refer  to  grid  indices  in  x-  and  y  directions  respectively, 
refer  to  partial  derivatives  with  respect  to  x  and  y  respectively 


refers  to  time  level  t 
refers  to  time  level  t+At 
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INTRODUCTION 

There  are  tvs  formulations  fur  the  numerical  solution  of  the  incompress ible 
Navicr  -  S  tokos  equations  ir  primitive  variables;  the  artificial  compress  ibi  1  i  ty 
ur.d  the  pressure  Poisson  methods.  In  both  methods,  the  velocity  field  is 
calculated  from  the  time  dependent  momentum  equation  using  time  march ing 
techniques,  while  each  method  employs  a  different  equation  to  compute  the 
pressure . 


In  the  artificial  compressibility  method,  s  time  derivative  of  the  pressure 
is  added  to  the  continuity  equation  [  1  ]  end  the  incompressible  field  is  treated 
es  compressible  during  the  transient  calculations.  On  the  ether  hand,  the 
pressure  Poisscn  method  ;2]  replaces  the  continuity  equation  with  a  second  order 
elliptic  Poisson  equation  for  the  pressure.  The  continuity  equation  it  then 
enforced  In  the  numerics. 


An  important  issue  in  the  numerical  solutions  of  the  primitive  variable 
formulations  is  the  satisfaction  of  the  discrete  continuity  equation  when 
staggered  or  non- staggered  computational  grids  are  used.  It  has  been  shown  for 
both  methods;  the  artificial  compressibility  ar.d  the  pressure  Poisson  that  the 
discrete  continuity  equation  is  satisfied  to  a  machine  zc:o  on  staggered  grids 
■  2,  3 J  .  Unfortunately,  this  is  net  true  on  non- stage tred  grids  The  artificial 
compressibility  method  requires  the  explicit  addition  of  a  fourth  order 
artificial  dissipation  term  to  the  discrete  continuity  equation  to  eliminate 
odd-even  decoupling  in  the  pressure  field.  The  odd-even  decoupling  is  caused  by 
central  second  erder  f inite -di f ference  approximations  of  the  first  order 
continuity  equation.  Therefore  the  discrete  divergence  of  the  velocity  field  is 
not  driven  to  a  machine  zero  but  rather  to  a  tens  proportional  to  the  fourth 
order  derivative  of  the  pressure  jdj.  Similarly,  the  pressure  Tcisson 
fcrmulation  does  not  satisfy  the  discrete  continuity  eouacion  exactly  on  non- 
staggered  grids.  This  phenomenon  can  be  explained  by  investigating  the  discrete 
pressure  equation  which  can  be  obtained  by: 


1.  Discretization  of  the  continuum  pressure  Poisson  equation  using  central 
second  order  accurate  formulas  f 5 .  6). 

2.  Direct  derivation  iron  the  discrete  divergence  of  the  discrete  momentum 
equation  ( 2 ] . 

Ai though  the  above  approaches  lead  to  the  sane  discrete  pressure  equation  on  a 
staggered  grid,  they  dc  not  give  the  same  equation  on  non-staggered  grids. 

Careful  examination  of  the  two  forms  of  the  discrete  pressure  equation 
show?  tnat  the  direct  derivation  of  the  discrete  pressure  equation  satisfies  the 
discrete  continuity  exactly  but  fails  to  give  smooth  pressure  field,  Th/e 
oscillatory  behavior  oi  the  pressure  is  caused  by  the  odd-even  decoupling 
inherented  in  the  resulting  discrete  pressure  equation.  Or  chu  other  hcr.d,  the 
first  approach  gives  smooth  pressure  field,  but  has  two  major  problems: 


The  compatibility  condition  of  the  Poisscn  Neumann  prcblen.  is  net 
automatical ly  satisfied. 


The  discrete  continuity  equation  is  not  exactly  satisfied- 


The  first  problem 
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The  equsti  ons  which  gc-'err,  the  laminar,  Ir.compress  ibie  flow  of  a  Nevtcr.ian 
fluid  .’re  given  in  Cartesian  coordinates  as  follows: 
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where  u  and  v  are  the  velocity  component:.  in  the  >;•  .od  y- directions 
rf5pH:;"i;y,  F  is  tbs  state  pressure  divided  by  the  density  and  ?.*.  is  the 
Rmvncids  -urine r . 

The  r.ai.n  difficulty  associated  with  the  solution  of  the  system  is 

the  continui ty  equation  fl).  Equation  (1)  is  a  constraint  which  the  velocity 
field  has  to  satisfy  at  any  instant  in  clrae  and  not  an  evolution  equation  of  the 
tyre  (2)  or  (3)  Also,  equation  (1)  does  not  involve  the  pressure,  which 
appears  only  in  the  nouentum  equations  (2)  and  (3). 

A  numerical  solution  for  equations  {!)•(})  must  incorporate  a  procedure  for 
taking  into  account  the  important  interaction  between  the  pressure  and  velocity 
fields.  In  order  to  achieve  this  coupling  between  the  pressure  and  tr.c  velocity 
fields,  the  continuity  equation  (1)  must  be  replaced  by  an  equation  which 
involves  both  pressure  and  velocity  and  at  the  seme  tire  guarantees  the 
satisfaction  of  the  incompressibility  constraint. 


.nq  Pressure  Poisson  rernu ration 


f;S  the 

the  momentum,  e 

the 

sh  .  ih 

ox2  5y2 

In  the  pressure  Poisson  formulation,  a  second -order  elliptic  eaustion,  cf 
Poisson  type,  for  the  pres-sure  is  derived  by  appjvfng  the  divergence  opt  rates  t>- 
the  momentum  equation  : 
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Che 


-,ui,  equation  r-fd-'ces  to: 
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(6) 


It  is  icportar.t  to  stress  here  that  the  same  re::, 
used  for  the  unsteady  reins  ir.  the  mc-men- turn  eq  :.?* 
i  nrc  rporn  ting,  equation  (6)  into  equation  we 


•r  di  y.r retisation  must  be 
.or.:.  ':)  and  (3).  By 


?  7 

C _ r  d _ F 
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Equation  (?)  is  the  pressure  Poisson  equation,  which  has  been  used  by  researchers 
to  resolve  incompressible  flows  on  staggered  ;  2  J  and  nor- staggered  grids  [6,5!. 
Tc  see  how  the  solution  of  equation  (7)  guarantee?  the  satisfaction  of  the 
continuity  equation  at  convergence-  wo  should  cerrsre  equation  (4)  with  equation 
(?)  rhis  comparison  reveals  that  the  pressure  Frisson  equation  (7)  can  be 
viev.-e.d  as  an  evolution  equation  for  the  dilation  D,  namely  equation  (6). 
Moreover.  Che  solution  of  equation  (6)  is  an  exponential  decay  in  time,  as 
follows : 


D(t)  -  D  e 
o 


■  t/fit 


(£) 


•-■here  D 


the  initial  dilation.  Thus,  the  solution  of  equation  (7; 


tor 


ores  sure .  along  with  equations  (2)  and  (3)  for  u  and  v,  guarantees  that  the 

initial  dilation  D  will  decay  to  aero  as  t—*  (i.e.  as  a  steady  state  is 
o 

approached) .  We  will  not  pursue  further  the  analvsis  of  the.  relation  between 
the  pressure  and  the  continuity  equations  in  continuum  form  •  a  subject  which 
has  been  investigated  by  Gresho  and  Sani  [7]. 

1 i sere t i ration  of  the  Governing  Equations: 

Using  the  Euler -explicit  temporal  di  sere 1 1  as t i on  scheme  for  the  tine 
derivatives  and  central  second  ordpr  accurate  firsts  difference  formulas  for  the 
spatial  derivatives.  th«  system  of  the  govern!,-.-  equations  (?)  ,  (7)  and  (3)  car¬ 
te  discretized  as  follows  (sc?  Fig.  1): 
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Vp  to  thin  point  cent:  nuity  equation  (1)  h->'  not  be'"  •  -  C  :■ 
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;c.f,  ti  r.'.'  i  ty  equation,  ’-bile  die  dilation  eZ  the  (t)  tire  level  is  retained, 
educes  tc: 
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important  to  stitts  here  that  the  same  temporal  disf.retitati  on  rust  be 


•jit.;  f.  r  the  or*,  teed-:  terms  in  the  momentur  equation?  (")  and  '1; 
;  roorp-e  rating  -.ousticr.  <;•• )  irt.o  equation  (o),  vc  obtain. 
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Equation  (7)  is  the  pre<suic  Poisscr.  equation  which  has  beer,  used  by  researchers 
to  resolve  ir.corprcsr.iblc  flow?  on  staggered  [2]  and  non- staggered  grids  [6,  S'). 
Tc  see  how  the  solution  of  equation  (7)  guarantees  the  satisfaction  of  the 
continuity  equation  at  convergence  we  should  cosoarc  equation  (l)  with  equation 
' 7)  .  Tl'.is  ccrrariton  rev-.-als  that  the  pressure  Poisson  equation  (7)  can  be 
vic-ved  as  ar,  ev.-iutior  equation  for  the  dilation  P,  namely  equation  (0). 
dereover,  the  solution  of  equation  (6)  is  an  exponential  decay  in  time,  as 
follows : 


D(t)  -  Do  e 


C/At 


(8) 


where  D  is  the  initial  dilation.  Thus,  the  solution  of  equation  (7)  for  the 
o 

treasure  ,  along  with  equations  (2)  and  (3)  fer  u  and  v.  guarantees  that  the 

initial  dilation  I>  will  decav  to  zero  as  t-»  (i  e.  as  a  steady  state  is 
o 

ipproached) .  Ve  will  not  pursue  further  the  analysis  of  the  relation  between 
the  pressure  ar.d  the  continuity  equations  in  continuum  ferr  -  a  subject  which 
her  been  investigated  by  Gresho  and  Sar.i  [ 7 j  . 

Ti sere t iza t ion  of  the  Governing  Ecuaflons: 

Using  the  Euler • explicit  temporal  discretization  scheme  for  the  tine 
derivatives  ar.d  centtf-i  second  order  accurate  finite  d. ffercr.ee  formulas  for  the 
spatial  derivatives,  the  system  cf  the  governing  equations  (7).  (2)  and  can 
se  discretized  as  fellows  (sc-c  Fig.  1): 
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Fig.  1.  Finite  Difference  Grid. 
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The  first  ar.d  second  order  ostraters  (5  ,6  )  and  (#  .  £  )  are  the  central 
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second  order  f )  ri te - difference  approximations  for  the  ,  )  and  ( — r.  * — r) 

cx  ty  ex2  <?■/ 

derivatives,  respectively.  The  pressure  equation  v.-.*  discretized  according  to 
the  consistent  finite  difference  ceded  proposed  by  Abdallah  [6]. 


ween 

-‘rich 


for  the 
3 '  can 


Equation  (?)  is  solved  for  ? 


sr.d  v.'  av.d  then  the  velocity 
'  j  } 


■  •  ,  given  u.  ... 

i.j  h)  1 

field  is  updated  using  equations  (10)  ar.d  (II).  At  this  point  of  our 
discussion,  it  seecs  important  to  pose  the  following  question:  Deep  the 
solution  c.f  equations  (9).  (10)  arc  fil)  guarantee  that  the  computed  velocity 
field  will  be  divergence  free  in  the  discrete  computational  space?  in  order  to 
answer  this  question,  ve  should  first  understand  how  the  discrete  oontinuitv 
equation  is  codclcd  ir  the  pressure  equation  (9).  Thf*.  point  car.  be  raade  very 
clear  if  wc-  look  at  the  derivation  of  equation  <9)  fret  a  different  point  cl. 
view,  hot  us  consider  the  following  discrete  approximation  cf  the  continuity 

to  satisfy  at  the  time  level  t+dt  ( n* : 


equation  (1;  which  we  seel 
n*l 
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The  pressure  Poisson  equation  '9)  can  be  derived  from  equation  (3  3)  If  vr  erpiov 
the  y  •  and  y-itoacrtur  "Ovation*.  at  the  nstfr  (iti/T.j)  and  (i 
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Subr.ti tur.ing  equations  ( j ua )  and  ''lib)  into  equation  (12;  wt  obtain: 
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It  takes  sort-e  simple  algebra  to  show  that,  equation  (15)  is  identical  to  equation 

'.?).  Notice,  that  the  tight  hand  tide  of  equation  (15)  is  the  divergent."  of  the 

velocity  field  while  its  left  hand  side  consists  of  the  x-  and  y-  components  cf 
the  steady  state  ferrr  cf  t:.e  merer,  turn  equation,  computed  at  the  nodes  (i±]/2.j) 
and  (i,j±l/2)  respectively.  Thus,  if  the  momentum  equations  are  driven  to 

steady  state  at  these  nodes,  the  left  hand  side  of  equation  (15)  eventually 

art-roaches  sere  and  rc  does  the  divergence  of  the  velocity  field 
unfortunate ly ,  cn  a  non- staggered  grid,  the  momentum  equations  are  not  driven  to 
tero  at  (i±l/2.j)  and  i'i.jnl/2)  nodes  but  they  arc  at  the  (i.j)  r.c-des  (see 
equations  (10)  and  (ID)  This  inconsistency  between  the  discrete  momentum  and 
continuity  equations  prevents  the  discrete  divergence  of  the  velocity  (right 
hand  side  of  equation  (15))  from  approaching  zzzo .  This  observation,  raises  the 
following  question:  what  is  the  sire  of  the  mass  source  (error  in  D)  which 
equation  (15)  introducer  into  the  flow  field?  Ve  answer  this  question  in  the 
following  section. 
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To  simplify  the  algibrc.  we  perform  the  following  analysis  using  the  Euler 
rather  than  the  Navier -Stokes  equations.  So,  we  drop  the  viscous  terms  from 
equations  (9),  (10)  and  (11).  This  simplification  does  not  alter  significantly 


the  generality  of  our  results 


For  the  soke  of  convenience,  let  us  introduce  the  following  notations: 
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Vich  reference  to  Fig  1.  we  obtain: 
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In  order  to  ®xprcsr.  ir  t  ttn-r.  of  the  d®pcr.d'.  nt  variable*  ’J .  %  a;-,'.’  p,  wC 
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"Die.  rig.*'--  hand  side  of  equation  (22)  is  the  rail  source  (error  in  the 
discrete  continuity)  which  is  introduced  into  the  i 2 o-  field  when  equation  (9) 
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Equation  !22)  reveal.*;  the  following  very  interesting  aspects  of  the 
pressure  Poisson  formulation  on  a  non-staggered  grid: 

i)  The  pressure  Polo son  method  satisfies  the  discrete  continuity  equation 
alcest  to  the  same  accuracy  as  the  artificial  compressibility  method  does. 
Recall  that  a  fourth  order  artificial  dissipation  term  is  explicitly  added 
to  the  pseudo-compressible  continuity  equation  of  the  artificial 
compressibility  method  ffi],  tc  stabilise  the  numerical  solution  In  the 
pressure  Poisson  formulation,  the  artificial  dissipation  is-  ’’implicitly" 
added  to  the  continuity  equation  because  of  the  way  the  pressure  equation 
is  discretircd  or.  a  non-staggered  grid. 

ii)  Cur  experience  with  the  pressure  Poisson  method  [6,  8]  and  ether 
researchers'  experience  with  the  artificial  compressibility  method  [a]  show 
that  the  dilation  can  be  quite  large  in  regions  of  the  flow  fiel.d  where 
high  pressure  gradients  occur.  High  pressure  gradients  result  in  high 
values  of  the  fourth  order  pressure  derivatives  .artificial  dissipation 
term.)  and  thus,  e  significant  mass  source  is  generated  in  the  continuity 
equation.  Fortunately ,  the  dilation  ir.  the  Foissrn  formulation  is 
explicitly  dependent  or.  the  square  of  the  grid  spiting  and  the  time 
increment.  Therefore,  to  control  the  errcrc  in  the  discrete  continuity,  wc 
strong's  recommend  th®  use  ef  fine  grid  ir  the  reri.nc  of  iigh  pressure, 
gradients . 
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Using  equations  (10,  and  (D)  into  ecustion  (24)  wc  obtain  the  following 
discrete  press  ire  equation: 
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•  -771  .  WR  r;n 
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(75b) 


P.y  inspecting  equation  (253).  we  car.  easily  conclude  that,  at  convergence ,  its 
left  hand  side  will  be  driven  to  zero  since  it  involves  the  steady  state  fern  of 
the  »  and  y-roruenttyn  equations  at  nodes  where  they  are  driven  to  zero. 
Therefore,  equation  (25;  will  satisfy  the  discrete  continuity  equation  to  a 
machine  tore. 


Unfortunately,  the  discrete  pressure  equation  (25)  produces  oscillatory 
solutions  for  the  pressure  because  of  odd-even  decoupling.  This  can  be  clearly 
seen  by  inspecting  the  discrete  Laplace  operator  in  the  left  hand  side  of 
ecustion  (25)  which  contains  either  odd  or  even  grid  points  in  the  x-  or 
y- directi  on.  Therefore,  in  rwo-dicensions  equation  (25)  gives  decoupled 
solutions  for  the  pressure  or.  the  odd-oed,  even-even,  odd  even  and  even-odd  grid 
points.  Each  solution  is  a  unique  solution  for  the  pressure  within  an  arbitrary 
constant.  Interestingly,  the  odd- odd  solution,  i or  example,  is  actually  a 
solution  for  the  pressure  eouaticn  cn  a  staggered  grid  which  is  a  sutset  of  the 
original  non- staggered  grid  (see  Fig.  2).  Similar  conclusions  can  be  drawn 
about  the  other  three  solutions. 
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In  conclusion,  the  use  of  equation  (9)  Co  compute  Che  pressure  will  produce 
a  smooth  pressure  find,  but  it  will  not  satisfy  the  discrete  continuity 
equation  exactly.  While  the  use  of  equation  (25)  will  produce  exactly  the 
opposite  result.  ;  r:  other  words,  on  a  non-stagger ed  grid  it  is  net  possible  t0 
se-tisfv  the  discrete  continuity  to  a  machine  zero  and.  at  the  same  tire,  obtain 
a  sro'-'-h  pressure  i:eld  T-.t  general.  one  r.as  to  sacrifice  partially  the 
satisfaction  cf  the  discrete  continuity  since  a  scocth  ar.d  physically  meaningful 
press  ire  field  is  desired. 

I r.  order  to  modify  equation  (25)  so  that  it  produces  snoot r>.  pressure  field, 
we  must  set  aside  the  idea  of  satisfying  the  discrete  continuity  to  machine 
accuracy  here  specifically,  we  seek  to  satisfy  the  discrete  continuity  up  tc  a 
fraction  of  the  mass  source  of  equation  (22).  as  follows 


n+1 


n+1 


v-l  . 
1  .  ’->1 


2  Ax  2/iy 

-  -  *  ii».?  ?  t?"  •) 

4  xxx;-:  I.j 


,n+l 

{  4-1 


-  ?-  i 
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yy>7 


(26) 


where-  t  is  a  positive  constant  ;0  <  <  <  1),  Usii-p,  the  rr.jncir.  -t  equations  (19) 
and  (11)  into  (25),  we  obtain  the  following  discrete  pressure  equation 


Lii-i^udL-iu  .  h.-r2r..t^ 


4*/ 
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ii±L_L 
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2  Ax 


± 

At 


J-LlLA  4  Iloiti 


(27) 


-*)••• 


22>v 


It  is  interesting  to  note  that  for  r  -  0  equation  (27)  reduces  to  equation  (?.ii) 
while  for  t-i  the  left  hand  side  of  equation  (27)  reduces  to  the  left  hand  side 
of  equation  (9)  furthermore,  both  equations  (l7)  and  (?)  p reduce  a  smooth 
pressure  field  but  they  fail  to  satisfy  the  discrete  continuity  to  ?  machine 
zero  However.  the  »rrot  in  the  discrete  continuity  which  is  prodyrtd  hv 
equation  '2')  is  r-ch  less  than  that  produced  by  equation  (9),  because: 

<J i  There  is  no  error  contribution  to  the  tars  source  (compare  equations 
(22)  and  (26))  from  the  source  term  of  the  pressure  equation  (27). 

(ii)  Numerical  experiments  show  that  the  odd-even  decoupling,  can  be  removed 
by  urir.g  numerical  values  for  <  «  l  (see  results  and  discussion 
section) . 


Anocher  advantage  that  equation  (27)  has  over  equation  (9)  is  associated  with 
the  computational  work  (CPU  tip®)  required  to  advance  the  pressure  and  velocity 
field®  t-  th®  rev  tire  level.  when  equation  (?'  is  solved  for  ?  along  with 
cat  -.cions  <19;  ar.d  (11)  for  u  and  v,  the  f  ano  r?  terns  need  to  be  discredited 
twice  ever--  tic-.-  step  at  ( i *  1/2 ,  jii/2)  for  the  pressure  equation  and  .?t 
(iil.ifi)  for  thi  moire n turn  equations.  This  is  not  the  case  when  equation  '"27) 
is  u--f*d,  rin'-c  b;‘h  the  pressure  and  mosier.tu.-B  equations  require-  the  calculation 
>f  the  f  tcror.  at  the  cam®  nodes.  Therefcre.  v<  re  cornice  nd  the  use  of 

:quwt{nn  27)  in  the  pressure  Poisson  formulations  or  non- staggered  gride,  since 
it  o.nim ;,7io  the  error  in  the  discrete-  continuity  (<  «  i)  and  reduces  the 
:or.' otior.al  tine  in  the  numerical  solutions. 
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Round-ary  Cor  die}  one-  for  th  e__  Pleasure  Equation: 

N'  uT^r.n  bounder/  conditions  for  the  pressure  Poisson  equarion  are  obtained 
>»sif.g  tr®  r.-rnul  component  of  the  aonentur  equation  along  the  boundary  contour. 


(28) 


condition  /28)  is  applied  at  tr.c  boundary  using,  one  sided  finite  difference 
•srrrexiir.atiens  for  the  x*deri‘.atives 
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(2?) 


••tcre  -  is  the  outward  unit  vector  normal  to  the  boundary  contour  S  enclosing 
the  solution  domain  A.  Equation  (9)  satisfies  identically  the  compatibility 
condition  on  a  r.on- staggered  grid  ( j  j  .  We  will  rhc>v  that  equation  (25) 
satisfies  the  integral  constraint  (2?)  as  well. 

For  the  sake  of  convenience  we  use  the  notation  introduced  in.  equations 
(It?)  and  0 6b) .  Incorporating  equations  (16)  into  (27)  we  obtain: 
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iiuiiliou  +  h^±± ;/huia  .  i  ,6x2  +  £i  26  j  ? 

2Ax  2 Ay  A  xxxx  y  yyy 


At 


(20) 


The  discrete  fern  of  the  integral  constraint  (29)  when  applied  to  equation  (30), 
reads  as  follows: 

br1  J;-5  .  hi.1n:hi.,i.i 

>2  j -i  ^  “>■ 


7  raxh  ..  *  4)-J, 

A  1  xx:-.  x  '  . 


„  IP;  . 

yyy>  *.j 


At 


-  0 


(31) 


••here  in  ar.d  jro  ere  the  numbers  of  the  maxi  suit  grid  points  in  the  x-  and 
y- directions  respectively.  Equation  (21)  can  be  simplified  as  follows. 

Af.  .  . -f,  <-f,  •  ia-I  h  .  J-h  .  .  , -h  . 
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(32) 


Arriving  ihf  :.>umanr  boundary  conditions  at  i-J  and  ir  (f,  ^  ^-0)  a"1^  stL 

>l  and  ;m  (h.  ^-h.  ^ 


,_-0)  and  noting  that  ,  j*hi  . ?”h- .  js- j 


. -0  at  stt^dv 


State  WC  can  eerily  see  that  the  first  two  terms  ir,  equation  (32)  vanish 
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i oent icfi - •* y .  FinA.ly,  tm.t  third  ter*  in  equation  (32)  is  als<  zero  sii.ee  the 
"continuity"  equation  is  locally  satisfied  at  every  node.  Theiefort,  equation 
(32)  and  consequently  equation  (21)  are  identically  satisfied  on  a  non -staggered 
grid. 


The  driven  cavity  7  r  cl  let  is  selected  here  as  the  model  problem  to  confir.T 
ou.-  analytical  <ir>  i  epm-er  is  .  Numerical  solutions  for  the  sr.oisentutr.  ^qu/stions 
(10)  and  (11)  are  obtained  rising  the  Euler  explicit  scheme  and  the  pressure 
I  tics  on  equation  using  'he  successive  over -relaxation  aethed.  Both  forms  of  the 
bis..:.: te  pressure  tquati'ns  <9,  and  (2~ }  are  vssd  to  invest  ig-ctc  their  ability 
to  satisfy  the  discrete  continuity. 

All  calculations  are  performed  on  uniform  non- staggered  grids  starting  from 

the  initial  conditions  u.  -  x.  .  .  v.  •  y.  and  P,  *»  0 .  The  reasc  r,  for 

1  j  i.j  i,j  1 ,  J  i.J 

this  choice  is  that  ve  vant  to  explore  the  capability  of  ti'.c  pressure  Poisson 
rethc-d  to  decay  an  initially  high  dilation. 

The  first  numeric?-  experiirent  is  conducted  ©"  a  ;  31x31)  grid  using; 
equation  as  the  discrete  pres.rure  equation.  The  artificial  dissipation 

parameter  <  is  set  to  zero.  The  time  evolution  of  the  logarithm  of  the  average 
dilation  at  a  r.cdc  is  shorn  in  Fig  3.  Cleat  ly  the  initially  high  dilation 
decays  to  the  machine  re.ro  (single  precision)  at  otcacy  state.  As  wc  expected, 
the  computed  velocity  f>e)<!  is  smooth,  vhiie  th--  pressure  ficid  is  oscillatory 
because  of  odd-even  decoupling. 


Convergence  of  D  (e  **  C) 
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I  ■  •  a  t  e  rj-  dec  ecol  lr.g  of  the  pressure  nodes,  ve  conducted  a  second 
the  sir..-  g;.:d  using  aquation  l''7)  with  non- zero  valuer-  c£  (.  Cur 
ate  that,  values  of  '.  as  low  .?s  0.1  are  sufficient  to  produce  smooth 
Luc.  The  logarithm  of  the  average  dilation  a:  *  point  is  shown  in 
-0.1  In  the  sane  figure  j.t  is  also  shown  the  dilattcr-  vhen 
is  used  to  compute  the  pressure  on  the  seme  grid.  As  car.  fca  seer 
:uarim  '27 '  w  th  £  -0.1  decays  the  logarithm  of  the  initial 
3  steady  state  value  0.005.  vh i ] e  the  corresponding,  value  for 
;s  0  02 .  The  results  confirm  our  analytical  developments  since  it 
that  equation  (9)  produces  h  hj.gher  dilation  than  equation  (27) 

1  We  should  mention  here  that  the  driven  cavity  problem  is 
difficult  case  35  far  as  the  satisfaction  of  the  discrete 
?  concerned.  The  reason  for  that  if.  the  existence  of  high  pressure 
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step,  (sec  Fig.  0)  is  used  as  u  second  test  case  in  order  to 
of  the  preceded  pressure  Pcrssc.n  solver  when  3nciyr.  irg  an 
h'.-  perforn  laminar  calculations  over  a  wj.de  range  of 
a  20ir.fl  uniform  grid.  The  Reynolds  number  is  bored  on 
t  «nd  the  downstream  boundary  is  located  at  a  -  1C h.  where  1 
The  computed  re  attachment  lengths  as  well  as  the  length  of 
tier  bubble  on  the  upper  vail  are  shown  ;n  Figs.  ~  and  £. 
omparisrn  with  the  eypetimentsi  data  of  F  c  f .  :9]  is  vc-ry 
For  hsgh'-r  Reynolds  numbers,  the  flow  field  becomes  three 
thus  cwo-Gimenr ional  calculations  are  not  sufficient.  The 
observed  in  other-  two  dissensions!  calculations  [10;,  jll) 
in  the  same  figures. 


U’c  showed  that  the  discrete  solutions  cf  incompressible  flows,  using  the 
classical  pressure  Frisson  formulation  cr.  non- staggered  grids,  do  ret  satisfy 
the  discrete  continuity  to  a  machine  accuracy  We  also  showed  that  the  error  i’ 
the  discrete  cent  ir.ui -y  i*  proportional  to  the  fourth  order  derivatives  of  the 
pressure,  the  trr.e  inc - oz?r,t  and  the  squere  of  the  grid  spacing.  It  is 
interesting  tc  point  cut  the t  the  ertor  is  proportional  to  the  explicitly  adoed 
artificial  dissipation  term  5.n  the  pseudo -cocpressibility  ms t hod  Sir.cc  the 
error  is  dt  i  e.vifn.t  upon  the  tine  and  spece  increments,  it  is  recommended  that  * 
fire  jr;d  ir  regions  of  high  pressure  gradients,  to  be  used  in  order  to 
ririmire  the  trier:-  It;  the  errors  ir.  the  discrete  continuity 
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Ve  derived  a  Oiodified  discrete  pressure  Foisson  equation  which  satisfies 
the  discrete  cortir.i  ty  cr  non-  staggered  grid-.  ur-  to  a  fraction  of  the 
dissipation  tens  Ir  the  classical  p  re  s  core  Foits.r  fortiulsci  on .  In  addition  to 
r ir miring  the  eirtr  In  the  discrete  continuity,  the  rrethod  also  reduces  the 
ci-cputational  wcrV  required  for  the  discretization  of  the  pressure  equation. 

f  Cr.h  G’-l.  EDG  Ft'  t>  7 

This  vc-rV.  was  spcrscred  by  the  Office  of  Naval  Research.  The  authors  thank 
:V  Von i tor  of  the  prc.crss  Mr  J  Fein,  for  the  support  cf  this  project. 
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